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A new approach to perform numerical simulations of systems with a ^-vacuum term is proposed, tested, and 
applied to CP 3 . The main new ingredient of this approach is the method used to compute the probability 
distribution function of the topological charge at 9 = 0. We comment also on the possibility of applying this 
approach to QCD at finite baryon density. 



The purpose of this contribution is to report 
a new approach, recently proposed by us in [Q, 
to simulate numerically systems with a ^-vacuum 
term. To fix the ideas let us consider the partition 
function Zy(ff) of any 9- vacuum like model in a 
finite space-time lattice volume V . This action 
can be written, up to a normalization constant, as 
the discrete Fourier transform of the probability 
distribution function (p.d.f.) of the topological 
charge at 6 = 0: 

Z V (9) = Y,Pv(n)e if>n , (1) 

n 

where pv (n) is the probability of the topological 
sector n. In almost all practical cases the sum in 
has a number of terms of order V since the 
maximum value of the topological charge at finite 
volume is of this order. 

The only a priori reliable numerical scheme 
to analyze the thermodynamics of ^-vacuum like 
models goes through the determination of the 
p.d.f. of the topological charge, pv(n), and the 
evaluation of (Q). But this is a difficult task due 
to the following two technical reasons: i) any nu- 
merical determination of pv (n) suffers from sta- 
tistical fluctuations Q , and small errors in pv (n) 
can induce enormous relative errors in the deter- 
mination of a quantity as Zy(&) which is an ex- 



tremely small number of order e~ v , ii) even if we 
were able to evaluate Pvip) with infinite preci- 
sion, the sum in (|l|) contains terms running from 
ltoe-^|. 

In few specific cases one can overcome the sign 
problem [Q . However previous attempts by other 
groups to simulate 9- vacuum like systems 
were based on the numerical determination of 
the p.d.f. of the topological charge straightfor- 
wardly, by standard simulations, or by more so- 
phisticated methods based on the use of multi- 
binning and re- weighting techniques. In all these 
attempts, artificial phase transitions at a 9 C de- 
creasing with the lattice volume were observed for 
the two-dimensional U(l) gauge theory at strong 
coupling as well as CP^ models. The origin of 
this artificial behavior, which follows from a flat- 
tening behavior of the free energy for ^-values 
larger than a certain 9 C , was analyzed in |^,^). 
Both groups agreed that the observed behavior 
was produced by the small statistical errors in 
the determination of the p.d.f. of the topological 
charge, the effect of which became more and more 
relevant as the lattice volume was increased. In 
Ref. H it was also noticed that by smoothing the 
p.d.f. flattening disappears. 

Our approach is based on a new method to 
compute the p.d.f. of the topological charge 
and the use of a multi-precision algorithm in or- 
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der to compute the sum in (|l|) with a preci- 
sion as high as desired. Let us write the par- 
tition function ([!]) as a sum over the density of 
topological charge x n = n/V and set py(n) = 
exp[— Vfv(x n )], where fv{x) is a smooth inter- 
polation of — 1/V \upv{x n ): 

Z v (8) = e- v ^^e WVx ". (2) 

We will assume that CP is realized in the vac- 
uum at 8 = since otherwise the theory would 
be ill-defined at 8 7^ @. This implies that 
exp[— Vfv (x n )] will approach a delta distribution 
centered at the origin in the infinite volume limit. 

Let us consider the partition function (Q) on 
the imaginary axis 8 = —ih, and let f{x) be the 
infinite volume limit of /y(x„). In the infinite 
volume limit the free energy is given by the saddle 
point. Assuming that f(x) has first derivative for 
any x except at most in isolated points, we can 
write the saddle point equation: 

f{x) = h. (3) 

Our proposal to compute the function f(x) is 
based in the following three steps: 

i. To perform standard numerical simulations 
of our system at imaginary 8 = —ih and to 
measure the mean value of the density of 
topological charge as a function of h with 
high accuracy (typically a fraction of per- 
cent). This is feasible since the system we 
have to simulate has a real action. Then, 
Eq. (||) is used to get a numerical evalua- 
tion of f{x). 

ii. To get f{x) we have to integrate f(x). Be- 
tween the possible ways to do this integral, 
we decided to fit f'(x) by the ratio of two 
polynomials, whose order is chosen to ob- 
tain a high quality fit, and then to perform 
analytically the integral of the fitting func- 
tion. In this way we get a very precise deter- 
mination of f(x), which allows us to com- 
pute the p.d.f. in a range varying several 
thousands orders of magnitude. This is the 
main advantage of our approach when com- 
pared with other methods based on a direct 
computation of pyin). 



iii. To use a multi-precision algorithm to com- 
pute the partition function (||) using as 
input the function f{x) previously deter- 
mined. 

The function f(x) obtained in step ii) suffers 
from statistical and systematic errors, the last 
coming from the fact that the saddle point equa- 
tion (|^) has finite volume corrections. An anal- 
ysis of these errors for the models and sizes we 
have studied, shows that systematic errors due to 
finite volume effects are smaller than the statisti- 
cal ones in the whole relevant range of x. This is 
the reason why we will replace fv(x n ) in equation 
(Q) by its asymptotic value f(x n ) in what follows. 
This substitution has no effect in the infinite vol- 
ume limit at imaginary 8 and we are assuming 
that the same holds for real 8. 

The x-interval used in the fits is an interest- 
ing point since it allows to improve our approach. 
From the definition of the partition function it is 
obvious that we do not need to know fv(x) for 
all the possible values of the topological charge 
density x: in the large volume limit only the val- 
ues of x such that fv{x) < gv{^) (where gv{8) is 
the exact free energy density normalized to zero 
at 8 = 0) are relevant for the computation of the 
free energy density. This allows in practice to 
fit only the relevant part of the data having very 
good fits with a suitably chosen low parameters 
function. 

To test these ideas we have analyzed three mod- 
els: the one-dimensional antiferromagnetic Ising 
model within an external imaginary magnetic 
field, the two-dimensional compact U(l) model 
with topological charge, and CP 3 in two dimen- 
sions. 

From a numerical point of view the determi- 
nation of the free energy density and order pa- 
rameter through equation (Q) in the Ising model 
has the same level of complexity of more so- 
phisticated models. Furthermore, in contrast to 
two-dimensional U(l) gauge theory, where the 
p.d.f. of the topological charge is nearly gaus- 
sian, the non-gaussian behavior of the p.d.f. of 
the mean magnetization in the antiferromagnetic 
Ising model makes this model a good laboratory 
to check the reliability of our approach. We refer 
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the interested reader to |l| for details of the sim- 
ulations, but would like to notice that the p.d.f. 
of the order parameter for such a system takes 
values in a range of around 2000 orders of mag- 
nitude. Notwithstanding that, we were able to 
reproduce the order parameter in the whole 6 in- 
terval within a few percent JlJ. The agreement 
between analytical and numerical results in com- 
pact U(l) was even better (within a few per thou- 
sand). Furthermore the flattening found in Q for 
the free energy density in relatively small lattices 
was absent in our simulations even in the 80x80 
lattice 0. 

The last model we have analyzed is CP 3 . It is 
the standard wisdom that this model shares many 
qualitative features with QCD. Even if it has not 
been analytically solved we believe it is worth- 
while to compare our results with previous ex- 
isting numerical simulations. Also in this model, 
the previous numerical simulations gave artificial 
phase transitions with a flattening behavior for 
the free energy density at a 9 C decreasing with 
the lattice volume (5jQ] . Our results for the order 
parameter at ft = 0.6 on a 100 2 lattice jl| have 
no trace of the fictitious phase transition found 
in (|] . Furthermore the order parameter is clearly 
different from zero at 6 = ir, hence CP is sponta- 
neously broken at this f3 value. An open question 
is how CP is realized in the continuum limit || . 

Finally we want to comment on the possi- 
bility of applying this approach to the analysis 
of QCD at finite baryon density. Indeed the 
partition function of QCD within an imaginary 
chemical potential belongs to the general class of 
partition functions described by (||), where x n 
now stands for the density of baryonic charge 
and exp[—V fv{x n )} is the probability distribu- 
tion function of the baryon density at vanishing 
chemical potential. The determination of this 
p.d.f. would allow us to solve QCD at finite 
baryon density. In this case we are involved with 
what we could call the inverse problem i.e., from 
numerical simulation of QCD at imaginary chem- 
ical potential one can measure the mean value 
of the density of baryonic charge and from it we 
can get the free energy density and the partition 
function of QCD at imaginary chemical poten- 
tial (left hand side of (§)). The last step of this 



procedure would be to compute, with the help of 
the multiprecision algorithm, the inverse Fourier 
transform of (j|) in order to get the p.d.f. of the 
density of baryonic charge. At the moment it is 
not obvious if this procedure would allow us to 
see the phase structure of QCD at finite chemical 
potential since, as follows from the discussion on 
the relevant x-interval in the determination of the 
topological charge as a function of 9, the range 
of density of baryonic charge at which the p.d.f. 
can be measured depends crucially on the maxi- 
mum value of the free energy density at imaginary 
chemical potential. However the field is so inter- 
esting that it is worthwhile to be investigated. 

Similar ideas to those presented in this contri- 
bution have been proposed and promisingly ap- 
plied to a random matrix model of finite density 
QCD in g. 
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